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Abstract 

We consider binary mixtures of soft repulsive spherical particles and calculate the depletion inter- 
action between two big spheres mediated by the fluid of small spheres, using different theoretical 
and simulation methods. The validity of the theoretical approach, a virial expansion in terms 
of the density of the small spheres, is checked against simulation results. Attention is given to 
the approach toward the hard-sphere limit, and to the effect of density and temperature on the 
strength of the depletion potential. Our results indicate, surprisingly, that even a modest degree 
of softness in the pair potential governing the direct interactions between the particles may lead to 
a significantly more attractive total effective potential for the big spheres than in the hard-sphere 
case. This might lead to significant differences in phase behavior, structure and dynamics of a 
binary mixture of soft repulsive spheres. In particular, a perturbative scheme is applied to predict 
the phase diagram of an effective system of big spheres interacting via depletion forces for a size 
ratio of small and big spheres of 0.2; this diagram includes the usual fluid-solid transition but, in 
the soft-sphere case, the metastable fluid-fluid transition, which is probably absent in hard-sphere 
mixtures, is close to being stable with respect to direct fluid-solid coexistence. Prom these results 
the interesting possibility arises that, for sufficiently soft repulsive particles, this phase transition 
could become stable. Possible implications for the phase behavior of real colloidal dispersions are 
discussed. 

PACS numbers: 82.70.Dd 
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I. INTRODUCTION 



Colloidal systems are complex fluids composed of mesoscopic particles dispersed in a 
solvent of microscopic particles. The large difference in size between the dispersed (colloidal) 
and the solvent entities may allow the latter to be considered as a continuun>i. 

One way to prevent colloidal particles from aggregating is to sterically stabilize them, i.e. 
to coat them with a suitable polymer layer. This leads to an effective repulsive interaction, 
frequently considered so steep that colloidal particles have been extensively modeled as hard 
bodies. Thus, a dispersion made of spherical colloidal particles may in a first approximation 
be assimilated to a system of hard spheres (HS), the phase behavior of which is well known. 
In this respect, systems of mono-sized polymethylmethacrylate particles provide a good 
example^. 

The case of binary colloidal dispersions, where two colloidal species of different size are 
immersed in a molecular solvent, is less clear. Here the HS model may also be used to char- 
acterize colloidal interactions. Sophisticated integral-equation theories for binary HS fluid 
mixtures indicate that, contrary to the prediction of the Percus-Yevick closure^, demixing 
between two phases of different compositions may occur, provided the diameters of the two 
components are sufficiently different^. Binary mixtures of silica particles were shown to 
undergo phase-separation^, but its nature was not ascertained. The possibility that the 
phases involved were both fluid could not be excluded, but recent work has given evidence 
that demixing involves fluid and solid phases, and that fluid-fluid segregation in HS mixtures 
should be excluded^*^'^'^. 

The mechanism supposed to be at the basis of these demixing phenomena is known as the 
depletion effect, according to which there exists an effective attractive interaction between 
the big particles due to their clustering leading to a larger free volume available for the 
small particles. It was first put forward by Asakura and Oosawa^*^, and later reconsidered 
by Vrij^^, to explain the phase behavior of colloid-polymer mixtures. Depletion effects are 
invoked to explain a large variety of phenomena in mixtures and have been the subject 
of intense theoretical effort in recent years. The validity of the effective one-component 
approximation (big colloidal particles interacting via a depletion potential) in describing the 
phase behavior of the actual binary HS mixture has been checked by comparison of the phase 
behavior obtained from the effective potential of mean force using computer simulation, with 
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that obtained from direct computer simulations of the mixture^i^. These studies revealed 
that fluid-fluid phase separation is actually unstable with respect to a fluid-solid phase 
transition; other computer simulation studies, based again on the effective one-component 
approach, showed no indication of any (even metastable) fluid-fluid demixing^^, suggesting 
that slightly different approximations for the depletion potential employed may lead to quite 
different phase behavior. 

The very shape of the depletion potential appears to be very sensitive to small variations 
in the direct interactions among particles. For instance, slight degrees of non-additivity 
in binary HS systems result in depletion potentials which are very different from those 
calculated for the corresponding additive binary mixtures--. Also, residual interactions of 
short range and moderate strength, such as very short-ranged attractive Yukawa tails, which 
prima facie would appear irrelevant compared to the dominant hard-body terms, have been 
shown to significantly affect phase behavior—, tending to make the effective interactions 
less attractive. Moreover, the importance of carefully taking into account the colloidal- 
colloidal direct interactions has been stressed^: the total effective potential is the sum of 
the direct and indirect (depletion) contributions, and clearly the final shape depends on very 
fine details of both. 

The above-mentioned sensitivity of (and on) the shape of the depletion potential un- 
dermines the rather common view that mixtures of sterically-stabilized spherical colloidal 
particles resemble HS mixtures, and serves as a stimulus for further investigations. In par- 
ticular, a largely ignored, and probably more realistic, variation on the hard-sphere theme 
is to consider particles to have soft, rather than hard, repulsive cores. Although quite 
convincing experimental evidence exists about the fact that colloidal particles coated with 
poly-12-hydroxystearic acid, one of the most frequently used coating agents, should really 
behave as hard spheres in many solvents^, it appears of relevance to investigate depletion 
interactions in repulsive particles of various degrees of softness. In fact, the softness of the 
interaction potential was shown to have a pronounced influence on the crystal nucleation of 
model colloidal systems^^, while interactions between colloidal particles may be tuned to a 
certain degree, for example by varying the chemical composition and thickness of the surface 
layer, to such an extent that soft-core colloidal particles can be realizec^^ . 

In this work, binary fluid mixtures of model soft repulsive spherical particles have been 
considered. The depletion interaction between two big particles induced by the small-particle 
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component has been evaluated using computer simulation. The depletion interaction can 
be obtained quite straightforwardly, using Molecular Dynamics (MD) simulation, as the po- 
tential of mean force on a big sphere. Systems of various densities have been analyzed and 
the results compared with various approximations based on a virial expansion. The most 
salient feature of the resulting depletion potential is that the effective attraction between 
soft spheres is greatly enhanced with respect to the HS the softness is increased, 

the attraction becomes substantially larger, and the effective diameter of the big spheres is 
reduced. Since these trends of the interaction potential tend to promote fiuid-fiuid separa- 
tion in liquids, the interesting possibility arises that fiuid-fiuid demixing could be stabilized 
in these fluids of repulsive particles. In an effort to explore this possibility, we have applied 
a perturbative scheme to evaluate the free energy, and from there the complete phase dia- 
gram, for a size ratio of small and big spheres of 0.2. These results already indicate that 
fiuid-fiuid demixing, which a related theoretical calculation predicts it to be absent in the 
case of HS mixtures, is close to being stable with respect to direct fiuid-solid phase separa- 
tion, suggesting that further variations of the interactions, for example by further softening 
the interaction potentials, might induce the stability of a liquid phase in mixtures of soft 
repulsive spheres. 

The article is organized as follows. In Sec. [TTl interaction models, simulation techniques 
and theoretical approximations for the depletion potential are described. Results are con- 
tained in Sec. IIIIl and the conclusions are presented in Sec. |Vl 

II. THEORETICAL APPROACH AND SIMULATION DETAILS 

The fluids considered in this work; are binary mixtures whose constituents are spherical 
particles of two species, i and j, interacting by a pair potential Uij{r) defined as 



with eij and aij the energy and range parameters, respectively, and r the distance between 
the two interacting particles, n, henceforth called 'softness parameter', is an exponent such 
that the limit n —>■ oo gives the HS case. The rate at which the HS limit is approached 
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as n ^ oo in Uij{r) is dictated by the temperature T. As n is decreased from this hmit 
the repulsive potential becomes softer. For the case n = 6 this potential corresponds to the 
repulsive part of the Weeks-Chandler-Andersen separation scheme for the Lennard- Jones 
potential^^. In this work, e.y=e Vi,j, while the size parameters for the two species have 
been given the following values: aii=a, cr22=0.2cr, (Ji2=cr2i=0.6a" (here an additivity rule is 
applied on the aij parameters; in this way the limit n ^ oo exactly corresponds to additive 
hard spheres). Particles of the species labeled as 1 will therefore be referred to as big, solute 
or colloidal spheres, while particles labeled as 2 will be called small or solvent spheres. 

One can define a total interaction potential U{r) between two isolated big spheres im- 
mersed in the fluid of small spheres. This potential will have two contributions: the direct, 
or bare, contribution, given by Mii(r), and an indirect contribution, W{r), which is trans- 
mitted through the small spheres and is actually the potential of mean force; this contri- 
bution is identified as the depletion potential, and the total potential will be calculated as 
f/(r) = Mii(r) + W{r). One of the aims of this work is to study how the depletion inter- 
action, and consequently the total effective interaction between the two big particles varies 
not only with the density of the solvent, p2, but also with the degree of softness, n, and the 
temperature T. 

In the case of generic particle systems, it is convenient to measure physical quantities in 
reduced units. Thus, we will define reduced densities as p* = Picrf^, i = 1,2, while a reduced 
temperature will be given as T* = k-BT/e = (/5e)^^, with j3 = l/k-QT and k-Q Boltzmann's 
constant. 

1. Virial expansion for the depletion potential 

The density functional formalism has been successfuly applied to calculate the depletion 
potential in HS mixtures^S"^. Due to the lack of density functionals of proven adequacy for 
soft-sphere systems^, we need to implement an alternative approach. One possibility is a 
systematic approach that provides a perturbative expansion for the depletion potential in 
powers of the solvent density. In the following we discuss such an expansion, which has been 
discussed before and applied to HS mixtures^i^^. 

Within the second-order virial approximation, the grand potential associated to the sol- 
vent particles, interacting via a pair potential U22{'r) and in the external potential Ui2{r) 
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corresponding to a colloidal particle located at the origin, is 

mp2] = J drp,{r) [lnp2(r) -l]+^-JdrJ dr' p,{r) p^ir') - r') 

+ P j c/r[Mi2(r) -p2]p2(r), (2) 

where P2(i') is the density distribution of the solvent particles, fij = 1 — exp {—(3uij) is the 
Mayer function associated with the interaction between particles of species i and j (note 
that here it is defined with a negative sign with respect to the standard definition). Within 
the same approximation, the chemical potential p2 corresponding to a reservoir of solvent 
particles with constant density p2 is given by 

/3p2 = lnp2 + P2 j dr'f22{r-r'). (3) 

Functional minimization of Eqn. ([2]) with respect to P2(r) and use of Eqn. ([3]) gives the 
following integral equation for the equilibrium density profile P2(r): 



P2[r) = P2exp 



-/?wi2(r) - J dr' [p2(r') - pa] /22(r - r') 



(4) 



We will obtain the solution to Eqn. (jlj) perturbatively, i.e. assuming the validity of the 
density expansion 

P2(r) =p2^i(r) + p^^2(r) + --- . (5) 

The functions \E'i(r) are found by inserting Eqn. into Eqn. (jl]); at first and second order 
in p2 we find 

^i(r) = exp[-/3Mi2(r)], (6) 
^2(r) = exp[-/3izi2(r)] j rfr7i2(r')/22(r - r'). (7) 

Up to third order in density, the contribution to the virial expansion of the excess part of 
the free-energy functional of the binary mixture containing terms linear in pi(r) is 

/3^W[pi,P2] = JdrJ dr'pi(r)p2(r')/i2(r-r') 

+ ljdrj dr' J rfr"pi(r)p2(r')p2(r")/i2(r - r')/i2(r - r")/22(r' - r"). (8) 

Using the definition of c^^'*(r; [pi, P2]), the one-body direct correlation function for big 
spheres, i.e. c{\r; [pi,P2]) = — 5/3J^ex[pi5 P2]/<^Pi(i"), we find the dilute limit pi ^ from 



Eqn. as 

-cS^V; [Pi = 0; P2]) = I dv'p,{v')h,{r - r') 

+ \ j dr' j dv"p^{r')p,{v")h^{v - r')/i2(r - r")/22(r' - r"). (9) 

The limit pi = represents the infinite dilution limit of the mixture with respect to the big 
spheres. As shown in Ref.— the depletion potential between two big particles in the sea 
of solvent particles can be obtained from 

PWiv) = cl') (oo; [pi = 0,p2]) - cf^(r; [pi = Cpa])- (10) 

Note that the calculation of the pair depletion potential requires the evaluation of the func- 
tional c^i\r] [pi = 0,p2]) at the equilibrium density profile of small particles in the presence 
of only one big sphere, P2(r), whose effect is taken care of by an external potential. This 
feature makes density functional theory to be a powerful tool to calculate solvation forces 
in multicomponent mixtures^^i^i. 

Substituting Eqn. ([9]) into Eqn. ( ITOl) and using Eqn. ([5]), with the expressions found for 
\E'i(r), Eqns. ([6]) and (171), we arrive finally at 

W{v) = iy(^o)(r) + W^(^=^)(r), (11) 
/?Vt^(AO)(r) = .p^l dr7i2(r')/i2(r-r'), (12) 

(5W^''^\v) = plj dr' j dr"Mr - r')/22(r' - r")/C(r, r', r"), (13) 
where the kernel 

/C(r, r', r") = e-"""''^''^ ^r") + ^ je-z^I^iaCrO+m^Cr")] _ f^^^^ _ r") (14) 

has been defined. The Asakura-Oosawa (AO) approximation given by Eqn. f[T^ . and the 
more accurate approximation given by Eqns. (|TT]) - (JT^ . which also includes a contribution 
from three-body interactions [see Eqn. ([8])] will be checked against computer simulations 
in Sec. IIIII All integrals involved in Eqns. ([T2l) and f|T3l) have been evaluated numerically 
using Gaussian quadratures. The above procedure can be systematically continued to obtain 
higher-order contributions to the depletion potential. However, inclusion of the next term, 
proportional to (p2)^, requires the evaluation of multidimensional integrals with a numerical 



8 



cost comparable to that required in MD simulations of the exact depletion potential. This 
is the reason why we stop the density expansion at second order in the present work. 

Let us now examine the predictions of the above two approximations with respect to 
density ^'^'^ softness parameter n. Fig. [T] shows, for values of the softness parameter 
n=6, 12 and oo, the depletion potential W{r) evaluated with the AO approximation, Eqn. 
( IT2|) . at a density = 0.0409, sufficiently low that the use of this approximation could 
in principle be justified, and at a temperature T*=0.776. The figure includes plots of the 
direct interaction Uii{r) between two bigger particles, as well as the resulting total interaction 
potential, f/(r) = uu{r) + W{r), obtained by the sum of direct and indirect contributions. 
Significant differences can be seen among the attractive depletion interactions as the softness 
parameter is changed: in particular, the softer the pair interactions, the longer-ranged the 
attractive contribution. As expected, the depletion potential in the case of the hard sphere 
binary mixture vanishes at a distance r* = r/cr = 1.2, while W{r) in the case of soft spheres 
does so at progressively larger distances. All depletion potentials follow the same trend, 
but the softest direct potential, corresponding to the case n = 6, exhibits a more attractive 
indirect interaction. This effect is partly counterbalanced by the fact that in this case the 
direct interaction potential is softer, remaining positive at larger distances than in the other 
cases (in fact, direct interactions vanish at a distance r* = 2^^/", which increases as n is 
reduced); the net effect is that the energy minimum f/min of the total interaction potential 
deepens and shifts to shorter distances as n is increased. Now, a convenient measure of 
the well width is given by the quantity Wn=(r^'^ — r^^)/r^\ with r^'^^ being the distance 
between the origin and the ith root of U{r). Using this criterion, it can be seen that the 
potential well becomes wider as the direct interactions operating in the mixture become 
softer. All these trends, which occur at the level of the AO approximation, are monotonic 
with the degree of softness. 

This situation changes to a certain extent when corrections from the second-order density 
terms are incorporated, as illustrated in Fig. [21 Here depletion potentials are calculated by 
employing Eqns. (III])-([T4]) at densities p^= 0.2045 [Fig|2](a)] and p^= 0.4090 [Fig. El^b)], 
both at T*=0.776. As expected, the depletion potentials become more negative with density. 
In addition, they start to show maxima at distances r* ~ 1.15-1.20; these maxima become 
more pronounced as density is increased. The existence of an associated repulsive barrier 
in the depletion potentials is a well-established fact, first shown in Ref.'^. The height of 
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the maximum is rather insensitive to n, while its location decreases very slightly with n. 
However, by looking jointly at Fig. 1 and 2, the combined action of P2 ^-^d n is noticeable. 
For relatively low values of p^-, the smallest value of U^i^ is observed for n oo, while 
f/inin is a monotonically decreasing function of n. In the same density regime, the width 
of the potential energy well, decreases steadily with p2 at a fixed value of n; the same 
monotonically decreasing behavior of Wn with n occurs even at higher densities, so that 
Wq > Wi2 > W(x, always. However, contrary to w^^,, which persists in decreasing with density, 
Wi2 at P2 — 0.4090 is almost the same as that at P2 — 0.2045, whereas wq at the higher 
density is even larger than the value at p^ = 0.04090. Finally, the value of Umin for finite 
n approaches that for n ^ 00, still the most negative, but in a nonmonotonic way, as the 
value of Umin ioT n = 6 is larger than that for n=12. 

Summing up to this point, we have shown that the density expansion of the depletion 
potential up to second order indicates that the predictions of the simple AO approximation, 
namely that the potential well associated to the total interaction shifts to longer distances 
and becomes less attractive with respect to the hard-sphere case as the direct interactions 
become softer, are different at higher fluid densities: the potential well and its depth be- 
have nonmonotonically with the softness parameter. A full calculation, valid to all orders in 
density, is needed to elucidate whether these trends could give rise to enhanced attractive in- 
teractions with respect to the hard-sphere mixture when the interactions become sufficiently 
soft. 

2. Details on MD simulations 

The predictions presented in the previous section are based on a density expansion of 
W{r). Since it is not practical to improve our approximation beyond the second order in 
density, one has to resort to computer simulation in order to check the validity of the results 
of the last section. We have carried out extensive MD simulations on a system of two big 
spheres immersed in a fluid of small spheres, calculating the depletion potential directly; 
this technique allows us to obtain depletion potentials which are, in principle, exact at all 
densities. 

In the above theoretical derivation, it has been assumed from the beginning that the two 
big spheres remained fixed at certain positions and that their effect on the small spheres 
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could be taken into account by means of an external potential. In our MD simulations, two 
methods have been used to freeze the dynamics of the big spheres. 

In the first method, the distance, r, between the big spheres has been maintained fixed in 
the course of the simulation run by using the method of constraints reviewed in Ref.— . Note 
that we assume the two spheres to form a dimer, with translational and rotational degrees 
of freedom; it is only the dynamics associated to the relative distance between the two big 
spheres which is frozen. This introduces an additional contribution to the free energy per 
particle with respect to the case where the centers of mass of the two big spheres are fixed, 
coming from the volume and solid angle explored by the dimer, but averaged quantities such 
as the potential of mean force are not affected. In the second method we fix the centers of 
mass of the big spheres, which therefore act as an external potential on the small spheres. In 
both methods, the force, F{r), acting on these bigger particles due to the smaller particles, 
has been calculated according to the formula first used in Ref.~ to evaluate the potential of 
mean force between two ions in a polar solvent, and already applied in the context of colloid 
physics in Ref.—, that is: 



In the above equation r is the unit vector defining the direction of the vector joining the two 



j. The angular brackets indicate an average over the configurations of the small particles 
generated during the simulation run. This method cannot be applied in HS mixtures, which 
require evaluation of the surface density distribution on the big spheres; possible inaccuracies 
introduced by extrapolation to contact of the density histogram obtained from the simulation 
are not present here. Statistical errors are of the same order, however, since in both, HS- 
and soft-sphere mixtures, the fluctuating quantity is averaged over 0{N2) values. 

Systems of two big particles, each of mass mi, immersed in a sea of A'"2=2046 or 4000 
small particles (depending on the method used to calculate the depletion force), each one 
of mass 7712 = have been simulated with standard MD simulation techniques. Periodic 
boundary conditions are applied on the rectangular simulation box. The simulations have 
been carried out for the cases n=6 and 12, integrating the equation of motion with the Verlet 
algorithn>2^, using a time step of 5xl0~'^r, r being equal to crii(mi/e)^''^, and maintaining 
the temperature at the desired value via the Nose-Hoover thermostat^. For each finite 




(15) 



big spheres, a and 6, while is the force on the big sphi 



.ere i = a,b due to the small sphere 
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value of n investigated, several distances have been considered, at which the depletion force 
is calculated during production runs of 1-3 x 10^ time steps, preceded by 5 x 10"^ time 
steps of equilibration. The depletion potential is then obtained by numerically integrating 
this force. Before integration, forces have been smoothed out by fitting to an appropriate 
function or by Akima spline approximation^-. The integration constant has been set so that 
the depletion potential was chosen to be zero at the largest radial distance considered. 

The solvent density, which should correspond to the density of the reservoir of solvent 
particles, is evaluated by computing the time-averaged density profile at the boundaries, 
i.e. far from the big spheres, averaging over the boundaries to suppress fiuctuations. This 
density does not exactly coincide with that obtained by dividing N2 by the volume available 
to the small spheres, since the density profile is highly structured in the region next to the 
big spheres. 

Using the techniques described above, depletion interaction potentials, W{r), and their 
corresponding total effective potentials, U{r), have been calculated for a few values of the 
density the softness parameter n and the temperature T*. 

III. PRESENTATION AND DISCUSSION OF THE RESULTS 

In this section we check the depletion potential resulting from the density expansion 
proposed above, see Eqns. (fTTll - (fT^ . against our simulation results. The purpose of this ex- 
ercise is to elucidate whether a truncated density expansion is adequate to describe depletion 
effects; in particular, we would like to estimate the maximum solvent density for which the 
truncated expansion is quantitatively valid. Then, the behavior of the depletion potentials 
obtained from simulation with respect to density, softness parameter and temperature will 
be assessed. 

A. Check of virial expansion against computer simulation 

First we compare the predictions obtained from the truncated virial expansion for a fixed 
softness parameter, n = 6, with the simulations results, at fixed temperature T*=1.000 and 
for several values of pg- Fig- [3] shows the results of this comparison. It appears that the first- 
order expansion, corresponding to the AO theory, insufficiently accounts for the simulation 
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data even at the lowest density considered, P2 =0.190. The incorporation of the next term 
in the density expansion improves significantly the performance of the theory with respect 
to simulation. The agreement between the outcome of Eqn. fllip and the simulation data 
is very good for =0.190 at all distances. Increasing the density, the agreement naturally 
deteriorates but is still quite fair at =0.381 for distances 0.80 < r* < 1.15. In contrast, the 
agreement between theory and simulation is poor at the very high density pg =0.740, which 
indicates that higher-than-second order terms should be included in the density expansion. 

Due to the "compensatory" effect of the addition of the direct interaction potential 'Uii(r) 
at the shortest distances, the agreement between the second-order virial expansion theory 
and simulation is better if the total effective potentials are examined (Fig. Hj). For p2 = 0.190, 
Eqn. ( ITTl) predicts a potential U (r) which is almost indistinguishable from that of simulation. 
At P2 = 0.381, the agreement between the two curves is still fairly good, especially with 
regard to the location and width of the attractive well, while the subsequent tiny maximum 
is displaced and its height slightly overestimated. It is only at larger values of the density of 
the smaller spherical particles, such as P2 = 0.565, that the two potentials show noticeable 
discrepancies. For pg = 0.740, the results disagree considerably even in the attractive region. 

By confronting Eqns. f|T2|) and ffTTHT^ together, and by exploiting the results just pre- 
sented, one can conclude that Asakura-Oosawa theory can be used to construct accurate 
total effective potentials up to density P2 — 0.05; the inclusion of the quadratic term is 
necessary for densities in the range 0.05 ^ P2 ^ 0.40, while for larger densities at least the 
cubic term is to be used. Once we have assessed the density range where the virial expan- 
sion is valid, we now turn to explore the dependence of the depletion and total effective 
potentials obtained from simulation with respect to density, the degree of softness and the 
temperature. 

B. Dependence of the depletion potential on density, n and temperature 

Fig. [5] shows the values of the depletion force, calculated from MD simulation, as a 
function of distance, for two values of the softness parameter, n = 6 and n = 12, and for 
two high densities, P2=0.680 and p2=0.750; the temperature is T* = 0.776. Note that the 
values of density chosen are certainly below the freezing density of the WCA fluid, which 
accurate estimations based on computer simulation locate at pg ~ 0.91 for T* = l.C'^''. In 
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all cases the depletion force has been evaluated in the region of small distances as well; 
note that F{r) vanishes at r = 0, i.e. when the two big spheres are exactly superimposed, 
because of symmetry. The computations of the depletion potential at such small distances 
could appear prima facie useless, as the direct interactions are so largely dominant in these 
conditions that the value of W{r) does not affect that of U{r). In particular, for hard-body 
systems, the evaluation of depletion forces at distances shorter than the contact distance is 
clearly useless in the construction of the effective pair potential. Indeed, previous simulation 
studies on hard-sphere binary mixtures have only reported data at distances larger than 
the contact distance^"^, and theoretical efforts have tried to reproduce these data {e.g. 
Ref.— 1^). Nonetheless, indirect interactions are clearly defined for all distances, regardless 
of the nature of the direct interactions, and their evaluation at very small distances should 
be of importance per se as they constitute a set of numerical data against which theories 
should be tested. 

The shapes of the depletion forces shown in Fig. [5] are clearly consistent with two distinct 
regimes as the radial distance r is changed: a shallow and deep minimum at scaled distances 
smaller than unity, and an oscillating part at longer distances. The value of the minimum 
is rather insensitive to the degree of softness but depends strongly on the value of density, 
while the location of the minimum is weakly n-dependent. Changing these two parameters 
gives rise to rather different oscillatory structures, the amplitude and damping of the latter 
being larger for increasing values of density and softness parameter. 

As shown in Fig. 5, the depletion forces have been fitted, using empirical functions 
which reflect the existence of these two regimes^. The depletion potentials have been 
then obtained by integrating the corresponding fitting function. The results are shown in 
Fig. O which contains the direct interaction potentials together with the resulting total 
effective pair potentials f/(r). Unfortunately, simulation data on the analogous hard-sphere 
binary systems are available only up to a density of P2 = 0.7^'^^. Thus, in order to make a 
comparison, the parametrized formula derived in Ref.— for the depletion potential in hard- 
sphere mixtures has been used at the highest density considered, = 0.750. From that 
figure, one can see that for the two densities considered, P2 = 0.680 and 0.750, not only the 
well width but also the depth of the total effective interactions between the big soft spheres 
are larger than in the hard-sphere case, and in the case of the softest potential, n = 6, the 
attractive interactions are considerably enhanced. 
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This result concurs with those of Ref.— about the importance of taking into account the 
direct interaction between the big particles. The total effective interactions which govern 
the thermodynamical, structural and dynamical properties of these particle systems, is in 
fact the sum of two contributions which may have, as in the present case, different signs but 
comparable absolute values in an interval of distances close to the range of the interaction 
between the big spheres. Thus, the final shape of the total effective pair potential is clearly 
very sensitive to the details of the direct and indirect contributions, both on the same footing. 

However, in contrast with the results in Ref.—, where attractive interactions were consid- 
ered and depletion interactions of varying nature (attractive or even repulsive) were obtained, 
in the case of repulsive direct interactions the depletion potentials are always attractive and, 
at high fluid densities, greatly enhanced with respect to those obtained for HS mixtures. In 
the light of these results, and considering that in some cases HS mixtures exhibit a metastable 
fiuid-fiuid phase separation, one may wonder whether enhanced attractive depletion forces 
may change the nature of demixing phenomena in mixtures. This topic will be discussed in 
the next section. 

One last issue in the present subsection involves the dependence of the depletion poten- 
tial on temperature, a dependence that certainly exists as the model particles interact via 
smooth, as opposed to hard, interactions. Therefore, the total interaction potential will in 
general be a parametric function not only of density, but also of temperature, U{r] P2,T*). 
Fig. 7 shows the MD-derived depletion and total effective potentials obtained for n=6 at 
densities P2 = 0.381 and 0.740, and temperatures T* = 0.776 and 1.000. As expected, 
density has a larger, dominant effect. Differences between the two depletion potentials 
corresponding to the same density but different temperatures are discernible only at the 
shortest distances, the two functions being essentially superimposed for r* > 1. In the case 
P2 = 0.381 the total interaction potentials are very close even at short distances, where the 
two functions are positive. When =0.740 the difference is larger at these distances, but 
still small. Irrespective of the value of density, the depletion potential becomes stronger as 
temperature is increased, which is reflected in the total effective potential having a deeper 
and wider attractive well; this can be understood in view of the higher momentum inter- 
changed between particles at higher temperatures, which involves larger kinetic pressures 
(and hence larger depletion effects) and more penetrable spheres. In this respect, although 
density can be considered in the present context as the more physically interesting variable, 
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being responsible for the largest variations in the effective potentials, the modest effect of 
temperature seen in Fig. [7] could have a large impact on the phase behavior and properties 
of binary mixtures of the type investigated in this work. 

IV. CALCULATION OF THE PHASE DIAGRAM FOR THE REPULSIVE SOFT 
SPHERE MIXTURE WITH n=6 

In the previous sections we have presented results for the depletion potentials of mixtures 
of soft spheres at various densities of the small spheres. These potentials have been com- 
puted, either using the virial expansion or by computer simulation, under the assumption of 
infinite dilution of the big spheres. Invoking the two-body approximation, we may assume 
that this is the potential felt by any two pairs of big spheres when their density is arbitrarily 
large. Indeed, the effect of three-particle interactions on the depletion potential has been 
argued to be overall very small^i^i^ii^^. In addition, we have made our own estimates of this 
effect by computing the force on a trimer of big soft spheres, in an equilateral triangular con- 
figuration, at a distance between the spheres r = an + a22, which approximately coincides 
with the maximum of the depletion potential, and then subtracting the pair force associated 
to the three bonds; the remainder is the three-body contribution. Our calculation indicates 
that this contribution can be neglected altogether, confirming previous expectations— i^>^i^ 
in hard-sphere mixtures. 

We now consider a collection of big spheres at density p* interacting via the depletion 
potential calculated at a density of small spheres pg- mentioned in the introduction, 
a longstanding discussion was whether or not a similar mixture of hard spheres exhibits 
demixing, i.e. a first-order phase transition where two phases with different concentration 
of big and small spheres coexist. The current understanding is that hard-sphere mixtures 
do exhibit demixing, but it involves a fluid phase and a crystalline phase, or two crystalline 
phases^. Since, as we have shown, soft repulsive spheres are influenced by a much stronger 
depletion force than in the case of hard spheres, it is pertinent to ask the following questions: 
(i) Does a mixture of soft-repulsive spheres exhibit demixing? (ii) Is the tendency toward 
segregation stronger than in the hard-sphere case? (iii) Does demixing involve two fluid 
phases, rather than one or the two being crystalline? Since we only intend to explore the 
possible phase behaviors of the mixture, we have applied a perturbation theory to calculate 
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the phase diagram; an alternative, in principle exact computer simulation approach is for 
the moment ruled out because of its relatively high computational cost. The theory is a 
density-functional theory based on the weighted-density approximation where a perturba- 
tion term is added, the so-called perturbative weighted density approximation (PWDA). 
The perturbative term is evaluated using the compressibility equation to fine-tune the ef- 
fective radial distribution function in the crystal phase. This procedure involves a mapping 
onto an effective HS fluid. Even though the theory cannot be trusted for liquid densities be- 
yond the crystallisation density for hard spheres, for densities less than that the predictions 
of the theory should be qualitatively correct. More details on the theoretical procedure can 
be found in Ref.— . 

The resulting phase diagram, in the p\-p2 plane, is shown in Fig. [HI The continuous 
lines correspond to the fluid-solid binodal lines obtained from the PWDA theory using 
depletion potentials calculated from the second-order virial expansion, Eqn. (fTTj) . The 
dashed lines are the corresponding calculation of the spinodal lines for fiuid-fiuid phase 
separation; note that the low-density spinodal should be given qualitatively correctly by 
the PWDA theory (maybe not so in the case of the high-density spinodal). It can be seen 
that fluid-fluid separation is always metastable with respect to the fluid-solid transition. 
PWDA calculations based on depletion potentials calculated from simulation are indicated 
by symbols. As expected, the two results are very similar for low values of where the 
truncated virial expansion is supposed to be accurate. At higher densities p2, the two sets 
of depletion potentials give different but compatible results (for both binodal and spinodal 
lines) . 

For the sake of illustration we include in Fig. [HI as dotted lines, the results for the 
fluid-solid binodal line from a slightly different version of perturbation theory^I as applied 
to the corresponding HS mixture using the depletion potential of Roth et al™. In fact, the 
HS fluid-fluid spinodal line predicted by this theory occurs at very high solvent densities, 
too high to appear in the range of Fig. [8] given that the effective depletion potential is 
not valid in this regime^, we conclude that, for the present value of size ratio, fluid-fluid 
phase separation in effective HS mixtures is, if anything, highly metastable with respect 
to direct fluid-solid phase equilibrium. In contrast, in the case of soft interactions, fluid- 
fluid phase separation is close to being stable with respect to fluid-solid phase separation. 
The predictions of perturbation theory as to the large differences between HS and soft- 
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sphere mixtures is very significant and illustrate the effect on the phase behavior of the 
significantly more attractive effective interactions operating in soft-sphere mixtures than in 
the corresponding HS mixtures. 

One can speculate on possible variations of the soft repulsive bare potentials (in the 
direction of making them even softer, for instance) that could promote the stability of the 
fluid phase, and its associated critical point, with respect to the solid phase. Actually, 
recent work— supports the possibility that an increasingly deeper and, especially, wider 
attractive well could lead to the observation of a stable fiuid-fiuid phase separation in these 
binary systems. In any case, use of an improved theory and more extensive simulations on 
fluids with other types of interactions and at different temperatures are needed, since the 
final outcome appears to be the result of a delicate free-energy balance between the fluid 
and the solid. Further exploration of these topics are left for future work. 

V. CONCLUDING REMARKS 

The present results are considered pertinent to real experimental binary colloidal disper- 
sions. It is known that many of them exhibit a fluid-fluid demixing transition and a phase 
behavior which is temperature dependent (for a discussion of the effect of temperature on 
the properties of sterically stabilized colloidal systems see Ref.~), two facts that are not cap- 
tured by assimilating the colloidal particles to hard bodies but that could be by modeling 
them as soft-core objects. 

Indeed, the pair-potential models investigated should be relevant to sterically stabilized 
colloids. In these systems, the degree of softness of the direct interactions can be controlled 
by varying the chemical composition and/or, especially, the thickness of the coating polymer 
layer: we could in principle expect that an increase in thickness produces softer interactions. 
Recently, an experimental study has been carried out on binary dispersions formed by col- 
loids sterically stabilized with a coating polymer layer of varying chemical composition and 
thicknes»i^. Although no stable fluid-fluid separation was observed, evidence for weakly 
metastable fluid-fluid demixing has been suggested for the mixture with the thickest sta- 
bilizing layer. These results have been interpreted in terms of a non-additive hard-sphere 
model. However, another, perhaps physically more sound, explanation of these experimental 
facts can be given in terms of the degree of interaction softness. An interesting sequel of 
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the study of Hennequin et al.— could be an investigation of additional binary suspensions of 
sterically stabilized colloids with varying thicknesses of the coating polymer layer. Further 
theoretical work along the lines presented in the previous section could also be valuable and 
guide possible experimental research avenues. 

We hope that the present results will stimulate the search for real colloidal binary mix- 
tures where interparticle interactions can be tuned from hard- to soft-sphere-like, as well 
as an investigation of how the phase behavior and the properties of such mixtures evolve 
upon changing such characteristics. Indeed, charge- and sterically stabilized monodisperse 
colloidal suspensions have recently been reported where the softness of the interactions is 
controlled by the concentration of the added salt^^. Although the pair potential models 
investigated in this work are not appropriate for charge-stabilized colloids, as the presence 
of charges requires a different functional form for the model pair interactions, the observed 
evolution of the depletion and total effective interactions with density and degree of softness 
of the direct interactions should be quite general. Thus, we believe it would be interesting 
to extend the study of Ref.— to bidisperse suspensions. 
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FIGURE CAPTIONS 



Figure 1: Dashed lines: depletion potentials calculated with Eqn. f|T2l) for n = Q, n = 12 
and n ^ oo (from bottom to top). Dotted lines: direct potential for n = Q (upper curve) 
and n = 12 (lower curve); the case n ^ oo is a vertical line at r = a. Continuous thick lines: 
from right to left, total effective potentials for n = Q, n = 12 and n oo. All calculations 
for a density P2 = 0.0409 and temperature r*=0.776. 

Figure 2: Dashed lines: depletion potentials calculated with Eqns. f|TTl) - f|T^ for n = 6, 
n = 12 and n oo. Dotted lines: direct potential for n = 6 (light gray) and n = 12 (gray); 
the case n ^ oo is a vertical line at r = cr. Continuous thick lines: total effective potentials, 
corresponding to softness parameters n = 6, 12 and oo indicated by labels and different gray 
intensity, (a) p* = 0.2045 and (b) p* = 0.4090. In both cases T*=0.776. 

Figure 3: Depletion potential for n=6 calculated using Eqn. f|T2l) (dotted lines), Eqns. 
(|TTl) -( fT4l) (dashed lines), and MD simulations (sohd lines) at densities P2 =0.190 (a), 0.381 
(b) and 0.740 (c). In all cases T* = 1.000. 

Figure 4: Total effective potentials for n=6 calculated using Eqns. f|TTl) - f|T4|) (dashed 
lines) and MD simulations (solid lines) at densities pg =0.190 (a), 0.381 (b), 0.565 (c) and 
0.740 (d). In all cases T* = 1.000. 

Figure 5: The depletion force as a function of distance for: n=6, P2=0.680 (a); n=6, 
p*=0.750 (b); n=12, p*=0.680 (c); n=12, p*=0.750 (d). Dots are the simulation results 
while lines are fitting curves^^. In all cases T* = 0.776. Insets are zooms of the regions 
about the maxima. 

Figure 6: Dashed lines: depletion potentials obtained by integrating the depletion force 
calculated in the simulations for n = 6 and 12 (value of n indicated as a label). Dotted lines: 
direct potentials for n = 6 and 12 (the case n = 6 corresponds to the softer potential). Solid 
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lines: total effective potentials for n = 6, 12 and oo (HS), indicated by labels and different 
gray intensities (light, medium and dark, respectively). Results in (a) are for p2 = 0.680, 
while those in (b) are for = 0.750. Data for HS mixtures (n = oo) from simulations of 
Ref.— (a), and from density functional theory of Ref.— (b). Insets are enlargements of the 
regions about the maxima. 

Figure 7: Effect of temperature on depletion (squares) and total effective (circles) in- 
teractions between two big spherical particles: (a) P2=0-381; (b) p2=0.740. Open symbols 
correspond to T* = 0.776, while filled circles were calculated with T* = 1.000. All results 
are for n=6. 

Figure 8: Phase diagram for the effective one-component system in the plane p\-p2- Solid 
lines: fluid-solid binodal lines, as predicted by PWDA theory using virial approximation for 
depletion potential, Eqns. f|TTl) - f|T^ : dashed lines: fiuid-fiuid spinodal lines as obtained 
from the same theory. Dotted lines: fluid-solid binodals for a hard-sphere mixture of the 
same size ratio, as obtained from perturbation theory (see Ref.—). Filled circles: fluid-solid 
phase boundaries as obtained from PWDA theory using depletion potential calculated from 
simulation. Open circles: fluid-fluid spinodal from the same theory. Filled squares: fluid- 
solid coexistence densities for the one-component WCA system from computer simulation 
(see Ref.—). The value of the reduced temperature is T* = 1.000. 
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